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Abstract 


A  new  semi-empirical  potential  in  iron  (Fe)  has  been  developed  based  on  the  quantum 
chemistry  concept-of  bond  order.  The  new  potential  was  calibrated  using  the  traditional  fitting 
to  the  universal  scaling  and  the  equilibrium  volume  and  cohesive  energy  of  base-centered  cubic 
(BCC)  iron.  With  a  total  of  15  fitted  parameters,  the  potential  reproduces  with  only  minor 
deviations  to  the  elastic  moduli,  the  volume-pressure  equation  of  states  in  the  BCC  phase,  the 
energies  in  face-centered  cubic  (FCC)  and  hexagonal  close  packing  (HCP)  modifications,  the 
BCC-HCP  phase  transformation  under  pressure,  and  the  energy  of  the  (111)  free  surface. 
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1.  Introduction 


Among  the  elemental  solids  which  have  received  special  attention  in  materials  theory  and 
applications,  iron  (Fe)  is  rather  unique  in  its  scientific  and  technological  importance.  These  areas 
of  importance  range  from  an  understanding  of  the  earth’s  inner  core  in  seismology,  to  structural 
applications  such  as  radiation  embrittlement  of  nuclear  reactor  pressure  vessels  (Wirth  et  al. 
1997)  and  toughness  of  high-strength  iron-based  alloys.  In  fact,  in  spite  of  emerging  families  of 
new  nonmetallic  materials,  iron  still  constitutes  the  technological  core  of  our  civilization. 

After  decades  of  research  on  iron  and  steel,  a  renewed  comprehensive  attack  on  this  material 
may  not  appear  particularly  timely  until  one  realizes  that  the  fundamental  understanding  of  many 
important  properties  (such  as  deformation  and  plasticity  of  base-centered  cubic  [BCC]  transition 
metals  or  chemisorption  of  gases  at  free  surfaces)  still  does  not  exist.  In  crystalline  solids,  the 
mechanisms  which  govern  the  bifurcation  of  behavior  between  brittle  and  ductile  responses,  as 
well  as  the  processes  at  free  surfaces,  reside  fundamentally  at  the  atomic  level. 

Ever  since  the  early  days  of  computational  modeling,  developing  the  capability  to  reliably 
predict  the  strength  and  other  physical  properties  of  iron  has  been  a  recognized  challenge.  In  the 
current  context  of  modeling  the  properties  and  behaviors  of  materials  across  various  length 
scales,  there  continues  to  be  interest  in  gaining  a  more  fundamental  and  quantitative 
understanding  of  how  iron  deforms  under  stress  and  behaves  under  other  physical  or  chemical 
exposures. 

Two  approaches  currently  exist  for  modeling  iron.  Development  in  the  recent  decade  of 
efficient  methods  of  ab-initio  calculations  and  proliferation  of  high-speed  computers  made  it 
possible  to  perform  high-precision  first-principles  calculations  on  iron.  This  approach  has 
essentially  improved  our  understanding  of  structural  properties  of  iron  under  high  pressure 
(Soderling  et  al.  1996),  as  well  as  elucidated  the  energetics  of  important  polymorphic 
transformations  in  iron:  BCC  face-centered  cubic  (FCC)  (Krasko  1987, 1989;  Krasko  and  Olson 
1989)  and  BCC-hardness  critical  process  (HCP)  (Eckman  et  al.  1998;  Sob  et  al.,  to  be  published). 


1 


Ab-initio  modeling  of  iron  grain  boundaries  has  resulted  in  an  understanding,  on  atom-electron 
level,  of  the  effect  of  interstitial  impurities  on  cohesion-decohesion  processes  at  grain  boundaries 
and  of  grain-boundary  embrittlement  (Krasko  and  Olson  1990;  Wu  et  al.  1994a,  1994b,  1996). 
The  state-of-the-art  ab-initio  calculations  on  iron-free  surfaces  has  enabled  one  to  better 
understand  the  multilayer  relaxation  and  chemisorption  of  hydrogen  (Wu  and  Freeman  1993a, 
1993b;  Krasko  et  al.,  to  be  published). 

However,  the  ab-initio  approach  at  the  present  time  has  two  fundamental  limitations.  Firstly, 
it  cannot  be  used  for  high  temperature  studies;  secondly,  modeling  processes  involving 
atomistically  large  systems  (even  of  a  few  hundred  atoms)  are  still  outside  of  our  computational 
capacity.  As  an  alternative,  certain  aspects  of  physical  behavior  can  be  simulated  at  the  atomistic 
level,  where  interatomic  interactions  are  described  by  empirical  (or  semi-empirical)  potential 
models.  However,  the  physical  rigor  of  such  results  may  be  uncertain  or  suspect  due  to 
inadequacy  of  the  potential  description. 

In  atomistic  simulations,  often  one  has  to  strike  a  balance  between  the  robustness  and 
accuracy  of  the  model  and  the  computational  feasibility  of  the  study.  BCC  transition  metals  pose 
particular  challenges  because  the  validity  of  using  empirical  classical  potential  functions  is  well 
known  as  questionable.  On  the  one  hand,  properties  such  as  strong  directional  bonding  and 
ferromagnetism  make  it  necessary  to  explicitly  consider  the  electronic  degrees  of  freedom  in 
describing  the  interatomic  interactions. 

Practically  all  current  atomistic  simulations  in  iron  are  concerned  with  BCC  structure  using 
variants  of  a  many-body  potential  in  one  of  two  forms;  the  so-called  embedded  atom  method 
(EAM)  (Baskes  et  al.  1988;  Daw  and  Baskes  1984;  Foils  et  al.  1986;  Daw  1989),  and  the  so- 
called  N-body  or  Finnis-Sinclair  (FS)  potential  (Finnis  and  Sinclair  1984,  1986). 

The  EAM  (Baskes  et  al.  1988;  Daw  and  Baskes  1984;  Foils  et  al.  1986;  Daw  1989) — the 
method  of  choice  in  virtually  all  atomistic  simulations — is  based  on  the  concept  of  substituting  a 
uniform  “embedding”  jelly  for  a  real  atomic  environment.  This  idea  originates  from  the  so  called 
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“effective  medium”  theory,  but  has  been  made  much  more  robust  by  incorporating  experimental 
information  for  fitting  the  model  parameters. 

According  to  EAM,  the  total  energy  of  a  condensed-matter  system  can  be  written  as 

E  =  Z,  Ffc>(R,)]  +  V2  2tJ  V^R,  -  Rj),  (1) 

where  F[p(R)]  is  the  “embedding  function,”  p(R)  is  a  quantity  characterizing  the  local  atomic 
environment  (e.g.,  the  free  atom  electron  charge  density  at  the  atomic  site  R  due  to  the 
surrounding  atoms),  and  Vrep(Ri  -  Rj)  is  a  repulsive  pair  potential.  In  the  original  EAM,  both  F 
and  V  were  model  functions  that  depended  on  a  number  of  parameters  being  fit  to  a  set  of 
experimental  data;  this  typically  includes  the  cohesive  energy  in  the  system  at  equilibrium,  the 
equilibrium  lattice  parameter,  the  elastic  moduli,  and  possibly  the  energies  of  metastable 
configurations  such  as  the  vacancy  formation  energy. 

Although  the  nature  and  shape  of  the  pair  potential  is  rather  clear,  the  electrostatic  repulsion 
(the  functional  form  of  the  embedding  function)  may  only  be  guessed.  Quite  a  few  different  forms 
of  F  and  other  ideas  have  been  suggested  since  the  EAM  was  formulated;  these  models  are 
sometimes  referred  to  as  “glue  models”  (Heine  and  Hafner  1990). 

As  mentioned  previously,  a  vast  family  of  EAM-type  models  for  iron  have  been  developed 
thus  far  (Simonelly  et  al.  1995).  The  most  advanced  models,  allowing  for  directionality  of 
interatomic  bonds,  actually  explore  two  ideas:  (1)  the  tight-binding  analysis  utilizing  moments  of 
the  electronic  density  of  states  (Carlsson  1991)  and  (2)  the  well-known  concept  of  bond  order 
(BO)  from  quantum  chemistry  (Horsfield  1996;  Horsfield  et  al.  1996;  Aoki  et  al.  1997;  Abel 
1985).  In  recent  years,  the  BO-type  empirical  potentials  have  been  widely  used  in  atomistic 
modeling  of  diamond-structure  semiconductors  and  their  alloys  and  compounds  (Tersoff  1988a, 
1988b,  1988c,  1989;  Tang  and  Yip  1995a,  1995b). 
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The  environmental  parameter  in  the  so-called  modified  embedded  atom  method  (MEAM) 
(Baskes  1997),  still  interpreted  as  an  electron  density,  has  a  nonspherical  symmetry,  which  allows 
for  directionality  of  bonds.  In  all  the  existing  EAM-type  models,  both  the  pair  potentials  and  the 
embedding  functions  are  of  a  more  or  less  arbitrary  character,  although  in  the  founding  papers 
(Daw  and  Baskes  1984;  Foils  et  al.  1986;  Daw  1989)  the  physical  connection  between  the 
embedding  function  and  the  exchange-correlation  energy  of  the  density  functional  theory  has  been 
established. 

Recently,  a  new  approach  explicitly  allowing  for  ferromagnetic  effects  via  the  Stoner  model 
has  been  suggested  (Krasko  1996;  Shibutani  et  al,  to  be  published).  In  this  report,  however,  we 
explore  the  ideas  behind  the  Tersofif  BO  approach  (Tersoff  1988a,  1988b,  1988c,  1989)  and 
develop  a  BO  potential  to  be  used  in  atomistic  simulations  in  iron. 

Section  2  briefly  outlines  the  BO  concept  and  the  Tersoff  potential  (Tersoff  1988a,  1988b, 
1988c,  1989),  which  served  as  a  prototype  for  our  new  potential.  In  section  3,  we  describe  this 
potential  and  the  procedure  of  fitting  its  parameters.  Section  4  is  a  conclusion. 

2.  Bond-Order  Concept  and  the  Tersoff  Potential 

The  BO  concept  follows  directly  from  quantum  chemistry.  Without  going  into  much  detail, 
we  will  briefly  outline  the  origin  of  the  BO  term  in  the  expression  for  the  cohesive  energy  (one 
can  find  a  concise  discussion  of  the  concept  in  Horsfield  [1996];  Horsfield  et  al.  [1996];  Aoki 
et  al.  [1997];  and  Abel  [1985]). 

From  the  tight-binding  or  linear  combination  atomic  orbital  (LCAO)  approach,  the 
self-consistent-field  Hamiltonian  is  written  as  (Abel  1985) 

H  =  T  +  IiV„  (2, 
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where  T  is  the  kinetic  energy  operator  and  V*  is  the  effective  potential  acting  on  atom  i.  It  is 
believed  to  be  short  ranged  and  strongly  localized.  The  summation  in  equation  2  is  extended  over 
all  atoms.  From  here,  an  expression  for  a  chemical  pseudopotential  immediately  follows,  with  the 
Schrodinger  equation 

{t  +  V,  +  X,  (v,  -  |<t>j  >  < <l>j|  V,)} |<>,  >  =  e,|4>,>.  (3) 

Here,  <f>i’s  are  the  atomic-like  local  orbitals,  and  V;  is  the  self-consistent  one-electron  potential  that 
includes  both  the  electrostatic  and  exchange-correlation  parts. 

From  this  starting  point,  the  expression  for  the  cohesive  energy  is 

£  =  1,(9,  Xj , ,  V„  (R„)  +  Xj , ,  b,V„  (Rs)}  (4) 

where  Vrep  (Ry)  is  the  repulsive  and  Vattr  (Ry)  is  the  attractive  interatomic  potential  between  atoms 
i  and  j  separated  by  distance  Ry.  These  potentials,  in  turn,  are  defined  as 

V„(Rlj)=«t,»i  |V]|fj>-«|,0i|4,“lx4,0,|V^“j>,and  (5a) 

v.»  (Rs)=  < <t>°i|Vi|^>0j  >,  (5b) 

where  are  the  “unperturbed”  atomic  basis  orbitals.  Factor  q;  is  the  electron  charge  density 
on  atom  i,  while  by  is  the  BO.  The  BO  is  the  strength  of  the  corresponding  (ij)  bond. 

BO  strongly  depends  on  environment.  The  number  of  neighbors,  Z,  is  the  simplest 
environmental  parameter.  It  was  shown  in  Abel  (1985)  that  at  large  Z  is 

b,  -  l/Vz.  (6) 
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Having  summed  over  the  nearest  neighbors,  j,  equation  4  will  result  in  the  attraction  term 

proportional  to  Vz,  which  is  exactly  what  the  FS  model  (Finnis  and  Sinclair  1984,  1986) 
postulates.  Equation  4,  however,  does  not  have  a  “glue”  form. 

In  Abel  (1985),  the  quantities  entering  equation  4  were  parametrized  by  simple  functions.  In 
the  recent  decade,  quite  a  few  expressions  for  the  both  environmental  parameter  and  the  BO 
function  were  suggested  (Balamante  at  el.  1992). 

Among  these  models,  the  so-called  Tersoff  potentials  (Tersoff  1988a,  1988b,  1988c,  1989) 
are  of  prime  interest  for  our  purposes.  Formulated  initially  for  silicon  crystals,  the  potentials 
were  also  used  for  carbon  and  germanium,  and  later  generalized  to  allow  for  multicomponent 
systems  description  (Tersoff  1988a,  1988b,  1988c,  1989;  Tang  and  Yip  1995a,  1995b).  The 
Tersoff  potentials  explicitly  take  into  account  the  directionality  of  bonds.  For  a  one-component 
covalent  crystal,  the  total  energy  is  taken  to  be 

E  =  1/2  E,„,  ft  (R,  -  R,)  [V„„  (R,  -  Rj)  +  VJR,^)] .  (7) 

Here,  fc(IRI)  is  the  cut-off  function,  where 

f.  <R|)  =  ■ 

rci  and  rc2  are  the  cut-off  radii.  The  repulsive  pair-wise  potential,  Vrep(Ri  -  R),  has  the  simple 
exponential  form 


1/2  +  1/2 cos  [jc(R  -  rcl )/( rc2  -  rcl )],  rcl  <  R  <  rc2 ,  and  (8) 

0.,  R  >  r.2 


Vrep  (R)  =  A  exp  (- p, |R| ). 


(9) 


The  adjustable  parameters  are  A  and  Pi.  Vb0(Ri,  Rj)  is  the  BO  potential,  defined  as 
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(10) 


vjR,.  R,)=  -B  b(R;,  Rj)  exp  (-p2|R,  -  R,|), 
where  b(Ri,  Rj)  is  the  “BO  function,” 


b(Ri,  Rj)  =  x[i  +  wr. 

(11) 

Cij  =  fc(|Ri  -  Rk|)g(0>  and 

(12) 

g(eijk )  =  1  +  c  I  -  l/|l  +  d  (h  -  cos6ijk  f  J  }. 

(13) 

Here,  is  an  environmental  parameter,  and  the  BO  function,  b(Ri,  Rj),  is  seen  to  behave  at  large 

Cas  l/VC- 


In  equations  9-13,  A,  Pi,  B,  p2,  P,  %,  n,  c,  d,  and  h  are  the  adjustable  parameters,  and  0iJfc  is 
the  angle  between  the  bonds  i  and  j  and  i  and  k. 

The  contribution  to  the  energy  of  the  attractive  BO  potential,  V,**  does  not  have  a  glue-model 
form.  However,  as  Brenner  (1989)  has  shown,  a  Tersoff  potential  with  special  parameter  values 
can  be  reduced  to  the  form  of  an  FS  model.  The  latter  is  not  surprising  since  the  glue-model  ideas 
directly  follow  from  the  quantum-mechanical  considerations,  as  does  the  BO  approach. 

We  will  describe  a  BO  potential  which  basically  utilizes  the  Tersoff  ideology  but  modifies 
TersofF  s  original  functions  to  make  the  method  applicable  to  transition  metals. 

3.  The  Bond-Order  Potential  for  Iron 

We  preserve  the  TersofFs  expression  for  the  total  energy  as 


E  =  y2  1,.;  ft(R,  -R,)|V„(R,  -  R,)  +  v^.R,)]. 


(14) 
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The  repulsion  part,  Vrep(Ri  -  Rj),  has  been  modified  to  allow  for  a  strong  repulsion  at  the 
distances  smaller  than  half  of  the  nearest  neighbor  distance  in  BCC  lattice,  R0,  where, 

Vrep(R)=  A  exp  (-Pj  |R|  +  P4(O.5R0/R)a).  (15) 

The  attractive,  BO  part  of  the  energy,  is  also  preserved  in 

vjR,,  R,)  =  -B[b(R,.  R,)]  exp  (-p|R,  -  R,|).  (16) 

However,  b(Ri,  Rj)  is  now 

b(R„  Rj)  =  C,/(l  +  Y,C„  +  YiCij2  +  Y,C„’  +  Y.Cii4  +  YsCij5  +  Y«C#‘).  (17) 

where 

Cij  =  ^k^i.j  fc(|Ri  -  Rk|)g(eijk)exp  (-p3|Rj  -  Rjl").  (18) 

We  use  the  same  cut-off  function,  (equation  10),  with  the  cut-off  radii  rci  and  rc2  between  the 
second  and  third  coordination  spheres,  as  in 


gfejj  =  1  +  c{l  -  p/[l  +  d(h2  -  coseijk2)2]}'5.  (19) 

Since  the  Tersoff  form  of  function  b(Ri,  R3)  was  quite  arbitrary,  with  the  only  requirement  of 
“correct”  (reciprocal  square  root)  behavior  at  large  £,  we  decided  to  try  a  more  general  form. 
The  exponent  f  was  allowed  to  take  both  positive  and  negative  values. 

The  adjustable  parameters  A,  pi,  B,  p2,  p3,  p,  c,  d,  8,  h,  y,,  and  y6nz  were  found  from  the 
following  conditions.  A  and  B  were  found  (at  fixed  values  of  the  other  parameters)  from  two 
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linear  equations  for  the  cohesion  energy  in  BCC  modification,  ECOh  at  the  experimental  equilibrium 
volume  Qo>  and  the  pressure  P  =  0  at  Q0  as  shown  in  Table  1.  These  equations  were  solved 
exactly  in  an  analytical  form. 


•  Table  1.  Experimental  and  Fitted  Properties  of  Iron  (Energies  in  eV,  Distances  in  A, 

Atomic  Volumes  in  A3,  Elastic  Moduli  in  102  GPa,  and  Surface  Energy  in 
eV/Surface  Atom) 


BCC 

Exp.  Value 

Calc.  Value 

a 

2.8589 

2.8589 

Qo 

11.6833 

11.6833 

-4.28 

-4.28 

||  Erep  =  2.11521;  Ebo  = -6.3952 

Cll 

2.4310 

2.3675 

1  C12 

1.3810 

1.3191 

C44 

1.2190 

1.2190 

C’ 

0.5250 

0.5242 

K 

1.7310 

1.6686 

Esurf(lll)a 

2.5 

2.2439 

||  FCC(£20  =  11.152)b 

||  Ecoh  -4.223 

-4.2229  || 

(Erep  =  2.4019;  Ebo  = -6.6249) 

HCP  (Qo=  10.398)c 

Ecoh  -4.213  -4.2134 

(Erep  =  2.9406;  Ebo= -7.1540) 

a  Estimated  based  on  ab-initio  calculations  (Wu  and  Freeman 
1993a,  1993b). 


b  The  atomic  volume  and  energies  as  found  from  experiments 
on  high-temperature  BCC-FCC  polymorphic  transformation 
(Bendick  and  Pepperhoff  1982). 
c  The  atomic  volume  and  energy  at  the  point  of  BCC-HCP 
phase  transformation  (130  kbar)  (Giles  et  al.  1971). 


Then,  the  multivariate  minimization  procedure  (IMSL  routine  DUMINF)  was  used.  We 
minimized  the  root  mean  square  (rms)  deviation  of  the  calculated  energy  from  the  following 
Smith-Banerjea  universal  scaling  law  (Banerjea  and  Smith  1988): 
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Eu,  =  Ecoh  0  +  s)exp  (-s), 


(20) 


where 


s  =  4.85406479  (-E„„/k  •  s0). 


(21) 


K  is  the  bulk  modulus  (in  Mbar),  and  So  is  the  Wigner-Seitz  radius  (in  A),  which  corresponds  to 
the  equilibrium  atomic  volume  Q0- 

Other  conditions  to  be  imposed  on  the  adjustable  parameters  were: 

•  The  calculated  elastic  moduli,  C\  and  C44  had  to  be  as  close  as  possible  to  their 
experimental  values  in  BCC  lattice.  The  bulk  modulus,  K,  would  be  automatically  close  to 
its  experimental  value  if  the  rms  universal  scaling  minimization  had  been  achieved. 

•  The  calculated  energies  of  FCC  and  HCP  modifications  at  the  corresponding  volumes  had  to 
be  as  close  as  possible  to  the  corresponding  experimental  energies  for  FCC  and  HCP 
lattices. 

•  The  calculated  (111)  free  surface  energy  for  BCC  lattice  had  to  be  as  close  as  possible  to  the 
corresponding  experimental  or  ab-initio  value. 

These  three  conditions  were  fulfilled  by  making  use  of  five  punishment  functions  equal  to  one 
if  the  experimental  and  calculated  quantities  to  be  fitted  were  exactly  equal;  the  function  would 
otherwise  be  growing  as  a  factor  times  the  corresponding  quadratic  deviation.  By  changing  the 
values  assigned  to  the  factors,  the  minimization  process  could  be  directed  towards  damping  down 
the  undesirable  deviation. 

Parameter  nz,  which  enters  the  expression  for  the  total  energy,  was  set  by  trial-and-error  tests 
nz  =  6.  We  “switched  off’  the  short  range  repulsion,  p4  =  0.  In  the  future,  a  more  appropriate 
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choice  of  the  value  for  p4  will  be  accepted;  however,  a  strong  increase  of  the  short-range 
repulsion  will  hardly  affect  the  energy  at  physically  reasonable  distances  but  may  be  important  in 
molecular  dynamics-type  atomistic  simulations  to  prevent  nonphysical  configurations  from 
occurring. 


Tables  1  and  2  present  the  results  of  the  parameter  optimization.  Table  1  shows  that  the  new 
potential  reproduces  the  elastic  moduli  with  very  high  precision.  It  is  well  known  that  C44  is 


Table  2.  Parameters  of  the  Potential 


A  =  0.2346 154E  +  04 

Pi  =  0.3465077E  +  01 

p4  =  0 

a  =  0 

B  =  0.1580257E  -  04 

P2  =  0.1 1 17197E  +  01 

b  =  Cf(l+YiC+.-Y6C6) 

f  =  0.4534376E  +  00 

Yi  =  0.3718635E  +  03 

Y2  =  0.1664427E  +  04 

Y3  =  -0.653982  IE  +  02 

Yt  =  0.9013512E  +  00 

Ys  =  -0.1142012E  -  01 

Y6  =  0.1385744E  -  03 

fc  (equation  5):  Cut-off  radii: 

rci  =  3.70,  rC2  =  3.60 

nz  =  6,  p3  =  0.6034363E  +  00 

g  =  1  +  c  {1  -  p/[l  +  d  (h2  -  cos62)2]  }(26) 

P  =  0.943367  IE  +  00 

5  =  0.9870393E  +  00 

c  =  0.419903E  +  01 

d  =  0.3752465E  +  02 

h  =  0.6915982E  +  00 
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especially  sensitive  to  directional  bonding  and  is  not  reproduced  well  by  the  traditional  EAM.  The 
energies  of  FCC  and  HCP  are  also  well  reproduced. 

An  interesting  feature  of  the  potential  is  that  the  exponent  f  in  the  BO  function  is  very  close  to 
0.5  (0.4534376),  rather  than  -0.5.  In  order  to  correctly  reproduce  the  energies  of  FCC  and  HCP 
modifications,  the  BO  function  has  to  increase  with  £  rather  than  decrease.  Our  attempts  to 
modify  the  function  so  that  it  behaved  correctly  at  large  £  failed. 

Figures  1-3  illustrate  some  of  the  results  of  the  parameter  fitting.  Figure  1  shows  the  total 
energy  as  a  function  of  the  volume  ratio  O/O0BCC,  (O0BCC  is  the  atomic  volume  of  BCC  iron  at 
equilibrium)  as  compared  to  the  universal  scaling  curve,  equation  14.  Figure  2  compares  the 
calculated  volume-pressure  equation  of  state  Q/£20BCC(P)  with  experimental  data;  the  agreement  is 
excellent. 

Fitting  the  energies  of  FCC  and  HCP  phases  of  iron  was  part  of  the  procedure  of  the  potential 
calibration.  The  HCP  experimental  energy  was  estimated  from  the  enthalpy  at  the  point  of 
BCC-HCP  phase  transformation  under  pressure  at  130  kbar.  Figure  3  shows  the  enthalpies  (H)  of 
BCC  and  HCP  phases  as  functions  of  pressure.  One  can  see  that  the  phase  transformation  occurs 
at  P  =  110  kbar.  This  should  be  considered  as  a  reasonable  agreement  with  the  experimental 
value  of  130  kbar. 


4.  Conclusion 

We  have  developed  a  semi-empirical  BO  potential  for  iron,  and  the  preliminary  results  look 
encouraging.  Although  some  further  testing  and  possibly  further  readjustment  of  the  potential 
parameters  may  be  needed,  hopefully  the  new  potential  will  be  instrumental  in  atomistic 
simulations  of  both  deformation  processes  in  iron  and  chemisorption  reactions  on  iron  surfaces. 


12 


Figure  1.  Comparison  of  the  Calculated  Volume  Dependence  of  the  Cohesive  Energy  With 
the  Universal  Scaling  (Eus)  Curve  (Banerjea  and  Smith  1998). 


0  0.05  0.1  0.15  0.2  0.25  0.3  0.35 

P(mbar) 


a  Giles  etal.  1971. 
b  Clendener  and  Drickamer  1964. 

Figure  2.  The  Calculated  Equation  of  States  for  BCC  Iron. 
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P(mbar) 

Figure  3.  The  BCC-HCP  Phase  Transformation  Under  Pressure. 

Recently,  a  method  of  explicitly  calculating  ferromagnetic  contributions  to  the  total  energy  in 
iron  based  on  the  Stoner  model  of  itinerant  ferromagnetism  was  suggested  (Krasko  1996; 
Shibutani  et  al.,  to  be  published).  If  successful,  this  method  may  become  an  ultimate  tool  for 
atomistic  simulations  in  iron.  While  the  detailed  implementation  and  testing  of  this  method  are 
still  in  the  future,  the  new  BO  potential  will  enable  valuable  information  on  some  important 
processes  in  iron  to  be  obtained.  A  BO  potential  of  the  same  type  will  also  be  used  in  the  Stoner 
model  method  for  evaluating  the  nonmagnetic  contributions  to  the  total  energy. 
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